#include <gmp.h>
#include <pthread.h>
#include <stdlib.h>
#include <stdio.h>
#define BUFFERSIZE 3
#define NUMITT 23

/*Este programa utilizará duas threads principais, A e Y.
  A thread A calculará os valores a<n> e a Y os valores de y<n>.
  A será dividida em mais duas sub-threads, A1 e A2 para calcular os
  dois termos da equação
*/

//Buffer que receberá os valores de Y<k>
typedef struct buff {

	mpf_t x[BUFFERSIZE];
	//mutex utilizado para a sincronização das threads.
	//Cada mutex bloqueia um espaço do buffer
	pthread_mutex_t mutex[BUFFERSIZE];
	//variável auxiliar
	mpf_t y0;
	//ponteiro do consumidor
	int ini;
	//ponteiro do produtor
	int fim;

} arg_buffer;

//Memória compartilhada entre as threads
typedef struct sh_mem {

	//termo a<n> do algoritmo de borwein
	mpf_t a;
	//buffer com histórico de termos y<n>
	arg_buffer y;
	//verifica se o y<0> já foi calculado
	//utilizado para executarmos a thread Y antes da thread A
	int begin;

} shared_mem;

//Dados para as subthreads
typedef struct data {

	//id da thread
	int sub_thread_id;
	//resultado dos calculos efetuados pela thread
	mpf_t result;
	//endereço da memória compartilhada
	shared_mem *sh_mem;	
	int k;


} sub_thread_data;

//inicializa o buffer de y<n> e configura seus valores iniciais
//arg_buffer *buff = endereço do buffer a ser inicializado
void init_buff(arg_buffer *buff) {

	int i = 0;
	//reseta o buffer
	for (i = 0;i < BUFFERSIZE; i++ ) {

		mpf_init(buff->x[i]);
		mpf_set_prec(buff->x[i], 33300000);
		mpf_set_d(buff->x[i],-1.0);
		pthread_mutex_init(&(buff->mutex[i]),NULL);

	}

	//configura valor de y<0>
	mpf_init(buff->y0);
	mpf_set_prec(buff->y0, 33300000);
	mpf_set_d(buff->y0,2.0);
	mpf_sqrt(buff->y0,buff->y0);
	mpf_sub_ui(buff->y0,buff->y0,1);

	//posiciona ponteiros no inicio do buffer
	buff->ini = 0;
	buff->fim = 0;

}arg_buffer *buff

//destrói o buffer liberando o espaço de memória ocupado pelo mesmo
void destroy_buff(arg_buffer *buff) {

	int i = 0;
	for (i = 0;i < BUFFERSIZE; i++ ) {

		mpf_clear(buff->x[i]);
		pthread_mutex_destroy(&(buff->mutex[i]));

	}

	mpf_clear(buff->y0);

}

//inicializa a memória compartilhada configurando seus valores iniciais
void init_shared_mem(shared_mem *sh_mem) {

	//y<0> ainda não calculado
	sh_mem->begin = 1;

	mpf_t aux;
	mpf_init(aux);
	mpf_set_prec(aux, 33300000);
	mpf_set_d(aux,2.0);

	mpf_sqrt(aux,aux);
	mpf_mul_ui(aux,aux,4);
	//aux = 4 * sqrt(2)


	mpf_init(sh_mem->a);
	mpf_set_prec(sh_mem->a,33300000);
	mpf_set_d(sh_mem->a, 6.0);
	//a = 6
	
	mpf_sub(sh_mem->a, sh_mem->a, aux);
	//a = 6 - 4 * sqrt(2)

	mpf_clear(aux);

	init_buff(&(sh_mem->y));

}

//destrói a memória compartilhada liberando seu respectivo espaço de memória
void destroy_shared_mem(shared_mem *sh_mem) {

	destroy_buff(&(sh_mem->y));
	mpf_clear(sh_mem->a);

	free(sh_mem);

}

//inicializa dados das threads de A (A1 ou A2)
void init_sub_thread_data(sub_thread_data *data, int id, shared_mem *sh) {

	data->sub_thread_id = id;
	mpf_init(data->result);
	mpf_set_prec(data->result,33300000);
	mpf_set_d(data->result,0.0);
	data->sh_mem = sh;

}

//libera espaço de memória utilizado pelas threads de A (A1 ou A2)
void clear_sub_thread_data(sub_thread_data *data) {

	mpf_clear(data->result);

}

//Calcula y<n>
void calc_y(shared_mem *sh_mem) {

	mpf_t aux, param, param2;
	mpf_init(aux);
	mpf_init(param);
	mpf_init(param2);
	mpf_set_prec(aux,33300000);
	mpf_set_prec(param,33300000);
	mpf_set_prec(param2,33300000);

	mpf_set_d(param,1.0); //param=1.0
	

	mpf_mul(aux,sh_mem->y.y0,sh_mem->y.y0); //aux = yk*yk = (yk^2)
	mpf_mul(aux,aux,aux);//aux = (yk^2)*(yk^2)

	mpf_sub(param,param,aux); //param = 1 - (yk^2)*(yk^2)
	mpf_sqrt(param,param); //param = sqrt(param)
	mpf_sqrt(param,param); //param = sqrt(sqrt(param))

	mpf_set_d(aux,1.0); //aux = 1

	mpf_sub(param2,aux,param); //param2 = 1 - sqrt(sqrt(1 - (yk^2)*(yk^2)))
	mpf_add(param,aux,param); //param = 1 + sqrt(sqrt(1 - (yk^2)*(yk^2)))

	mpf_div(sh_mem->y.x[sh_mem->y.fim], param2, param);
	mpf_set(sh_mem->y.y0, sh_mem->y.x[sh_mem->y.fim]);

	int prev = sh_mem->y.fim;
	// Produz 1 elemento no buffer, posiciona o ponteiro de fim para o proximo elemento
	sh_mem->y.fim = (sh_mem->y.fim + 1) % BUFFERSIZE;
	
	//Lock na RC posterior
	pthread_mutex_lock(&(sh_mem->y.mutex[sh_mem->y.fim]));
	//Unlock na RC atual
	pthread_mutex_unlock(&(sh_mem->y.mutex[prev]));

	//inicializa a execução da thread A
	sh_mem->begin = 0;
	

	mpf_clear(aux);
	mpf_clear(param);
	mpf_clear(param2);
	
}

void *monitor_y(void *args) { //producer

	int i = 0;

	shared_mem *sh_mem = (shared_mem*) args;
	
	//Lock na RC 0
	pthread_mutex_lock(&(sh_mem->y.mutex[sh_mem->y.fim]));

	for ( i = 0; i < NUMITT ; i++ ){

		calc_y(sh_mem);
			
	}

	//libera ultima RC para thread A poder ser finalizada
	pthread_mutex_unlock(&(sh_mem->y.mutex[sh_mem->y.fim]));
	pthread_exit(NULL);

}

//Calcula o primeiro termo de a<n> -> a<k> * (1 + y<k+1>)² * (1 + y<k+1>)²
void *a_first_term(void *args) {

	mpf_t param;
	mpf_init(param);
	mpf_set_prec(param,33300000);

	sub_thread_data *sb_data = (sub_thread_data*) args;

	mpf_add_ui(param, sb_data->sh_mem->y.x[sb_data->sh_mem->y.ini], 1.0); //param = 1 + y<k+1>
	mpf_mul(param, param, param); //param = (1 + y<k+1>)²
	mpf_mul(param, param, param); //param = (1 + y<k+1>)² * (1 + y<k+1>)²
	mpf_mul(sb_data->result, sb_data->sh_mem->a, param); //result = a<k> * (1 + y<k+1>)² * (1 + y<k+1>)²

	mpf_clear(param);

	pthread_exit(NULL);

}

//Calcula segundo termo de a<n> -> 2 ^ (2 * k + 3) * y<k+1> * (y<k+1> + 1 + y<k+1> * y<k+1>)
void *a_second_term(void *args) {

	mpf_t param, param2, y2;
	mpf_init(param);
	mpf_set_prec(param,33300000);
	mpf_init(param2);
	mpf_set_prec(param2,33300000);
	mpf_init(y2);
	mpf_set_prec(y2,33300000);

	mpf_set_d(param, 2.0);
	mpf_set_d(param2, 1.0);
	mpf_set_d(y2, 0.0);

	sub_thread_data *sb_data = (sub_thread_data*) args;
	
	mpf_pow_ui(param, param, ((2 * sb_data->k)+3)); //param = 2 ^ (2 * k + 3)

	mpf_mul(param, param, sb_data->sh_mem->y.x[sb_data->sh_mem->y.ini]); // 2 ^ (2 * k + 3) * y<k+1>
	mpf_add(param2, sb_data->sh_mem->y.x[sb_data->sh_mem->y.ini], param2); //param2 = y<k+1> + 1
	mpf_mul(y2, sb_data->sh_mem->y.x[sb_data->sh_mem->y.ini], sb_data->sh_mem->y.x[sb_data->sh_mem->y.ini]); //y2 = y<k+1> * y<k+1>
	mpf_add(param2, param2, y2); //param2 = y<k+1> + 1 + y<k+1> * y<k+1>
	mpf_mul(sb_data->result, param, param2); 
	//result = 2 ^ (2 * k + 3) * y<k+1> * (y<k+1> + 1 + y<k+1> * y<k+1>)

	mpf_clear(y2);
	mpf_clear(param);
	mpf_clear(param2);

	pthread_exit(NULL);	

}//Lock na RC posterior

//Calcula a<n>
void calc_a(shared_mem *sh_mem, int i) {

	/*
		st_data[0] = ak * ((1 + y<k+1>)^4)
		st_data[1] = (2^(2 * k + 3)) * y<k+1> * (1 + y<k+1> + y²<k+1>) 
	*/
	sub_thread_data st_data[2]; 
	pthread_t threads[2];

	//inicializa os dados das threads
	init_sub_thread_data(&st_data[0],0,sh_mem);
	init_sub_thread_data(&st_data[1],1,sh_mem);

	//iteração de a<k>
	st_data[0].k = i;
	st_data[1].k = i;

	pthread_create(&threads[0], NULL, a_first_term, &st_data[0]);
	pthread_create(&threads[1], NULL, a_second_term, &st_data[1]);

	pthread_join(threads[0],NULL);
	pthread_join(threads[1],NULL);

	mpf_sub(sh_mem->a, st_data[0].result, st_data[1].result); 
	//a<k+1> = a<k> * (1 + y<k+1>)² * (1 + y<k+1>)² - 2 ^ (2 * k + 3) * y<k+1> * (y<k+1> + 1 + y<k+1> * y<k+1>)

	// Consome 1 elemento do buffer, posiciona o ponteiro de inicio para o proximo elemento
	int prev = sh_mem->y.ini;
	sh_mem->y.ini = (sh_mem->y.ini + 1) % BUFFERSIZE;
	
	//Lock na RC posterior
	pthread_mutex_lock(&(sh_mem->y.mutex[sh_mem->y.ini]));
	//Unlock na RC atual
	pthread_mutex_unlock(&(sh_mem->y.mutex[prev]));

	clear_sub_thread_data(&st_data[0]);
	clear_sub_thread_data(&st_data[1]);

}

void *monitor_a(void *args) { //consumer

	int i = 0;

	shared_mem *sh_mem = (shared_mem*) args;
	
	//enquanto a<0> não for calculado, espere.
	while(sh_mem->begin);

	//Lock na RC 0
	pthread_mutex_lock(&(sh_mem->y.mutex[sh_mem->y.ini]));

	for ( i = 0; i < NUMITT; i++ ) {
		calc_a(sh_mem,i);
	}

	pthread_exit(NULL);

}

int main() {

	//memória compartilhada
	FILE *f;
	int rc;
	shared_mem *sh_mem = (shared_mem*)malloc(sizeof(shared_mem));
	init_shared_mem(sh_mem);
	pthread_t threads[2];

	mpf_t pi;
	mpf_init(pi);
	mpf_set_prec(pi, 33300000);
	mpf_set_d(pi,1.0);

	//Cria Thread Y
	rc = pthread_create(&threads[0], NULL, monitor_y, sh_mem);
	
	if (rc) {
		//printf("ERRO AO CRIAR A PRIMEIRA THREAD\n");
		exit(1);
	}

	//Cria Thread A
	rc = pthread_create(&threads[1], NULL, monitor_a, sh_mem);

	if (rc) {
		//printf("ERRO AO CRIAR A SEGUNDA THREAD\n");
		exit(1);
	}

	//Espera Threads A e Y
	pthread_join(threads[0], NULL);
	pthread_join(threads[1], NULL);

	mpf_div(pi,pi,sh_mem->a);
	
	//gmp_printf("%F1.50f\n",pi);
	
	f=fopen("saida.txt","w");
	gmp_fprintf (f, "%F1.10000000f\n",pi);

	destroy_shared_mem(sh_mem);

	return 0;

}


